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Abstract This lecture describes the treatment of heavy 
quarks in lattice QCD by implementing the Isgur-Wise limit. 
The method is briefly discussed, and some of the special 
features of the resulting theory are highlighted. We 
emphasize issues of the renormalization of the effective 
theory. The formulation permits a calculation of heavy 
quark processes even when the momentum transfers are 
much larger than the inverse lattice spacing. Applications 
include semi-leptonic heavy quark decay and scattering 
processes, including the computation of the nonperturbative 
part of the Isgur-Wise universal function. 



1. INTRODUCTION 

Heavy quarks are an especially interesting challenge in lattice 
QCD. In many ways heavy quarks should be much easier to treat than 
light ones, since their Fermi momentum is typically non-relativistic. 
However, the three-momentum of even a slow-moving heavy quark can 
easily be larger than 2n/a (where a is the lattice spacing), the maximum 
value of each momentum component that is representable on a lattice. 
A first glance, then, the whole idea of discussing heavy quarks in a 
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lattice context seems problematical. Fortunately, for some problems we 
can insure that the quark's momentum is small enough so that it is a 
good approximation to take a strictly static heavy quark as a starting 
point, and develop an expansion in either its inverse mass, its velocity, 
or its momentum. The striking success in calculating the spectrum of T 
resonances, which are approximately bb bound states, shows the power 
of this approach when the momenta are under good control. On the 
other hand, there are important processes in which one cannot constrain 
the momenta to be small. 

A basic process which necessarily involves dealing with large 
momenta is the decay of one particle containing a heavy quark into 
another, such as the semileptonic B meson decay, B * - Dlv . In the 
continuum there has been great progress in efficiently describing the 
dynamics of such processes. There is a sensible expansion in the inverse 
of the heavy quark mass, and for infinite mass, known as the Isgur-Wise 
limit, notable simplifications obtain. The analysis to be presented here 
concerns the implementation of this dynamical limit on the lattice. We 
begin by briefly reviewing the continuum heavy quark theory, and the 
lattice theory which describes it. We then discuss the renormalization of 
the resulting Euclidean-space theory, and introduce a method for 
efficiently performing lattice calculations with it. We conclude with 
short discussion of some of the more interesting problems that are 
probably ripe for solution. 



2. THE ISGUR-WISE LIMIT 

The key dynamical idea of the Isgur-Wise limit is that the four- 
momentum of a heavy quark can be usefully decomposed as 
= Mv^ +p )i , where v is a constant but arbitrary four-velocity 
normalized to v 2 = 1 [1][2][3]. The four-vector v is called the classical 
velocity of the heavy quark. In the M - °° limit, a new symmetry 
emerges. Part of the symmetry is a flavor SU(N) symmetry between N 
flavors of heavy quarks. It relates, for example, decays of B mesons to 
those of D mesons. The more powerful aspect of the symmetry comes 
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from the fact that in the M - °° limit, the spin of the heavy quark 
becomes dynamically irrelevant. This aspect of the new symmetry 
relates processes involving particles of different spin, such as B - D 
and B - D* transitions. 

The simplifications of the Isgur-Wise limit are most directly 
seen in the perturbative fermion propagator. Assume we have a heavy 
quark of mass M, moving with momentum P. If we write P as Mv + p 
and pass to the M - °° limit, the fermion propagator simplifies to 



1 1+yv 

(1) 



y-P-M 2vp 

In the limit the heavy quark mass M and the large momentum Mv 
disappear, leaving behind only finite terms. 

The equivalent to the Isgur-Wise limit in Euclidean coordinate 
space is the following. The propagator is given by 



S( X ) = | J^P- _L_ e*' (2) 
1 (2ti) 4 iyp+M 



which satisfies the equation 

iyd+M)Six) = 8(x) (3) 

By an appropriate shift of contour in the energy plane, one can show that 
in the M - °° limit, for positive Euclidean time, 



e M VA - iM ^ ^ d*p_ 1-iyv eip . x 

' (2ti) 4 2p-v {V 
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where v is the four- vector (iv Q ,v) and v Q is [l+v 2 ] 1/2 . For x 4 < 0, 
the limiting form of the propagator is the same but with the sign of v 
(and v 4 ) reversed. This represents the heart of the Isgur-Wise limit 
because there are no kinematically large momenta remaining in the 
integral. The appearance of the complex energy -momentum four-vector 
is unusual. It ensures that the propagation of a heavy quark moving with 
three-momentum Mv shows a decay in Euclidean time with the correct 
energy, and that the projection operator appearing in the numerator has 
the same form it has in Minkowski space. 

The heavy quark limit can be compactly expressed in terms of 
an effective Lagrangian. One defines a reduced field by factoring from 
the quark field the phase which diverges in the M -> °° limit, 

h (v) (x) = lim e~ iMvx 1 + Y ' V i|/(x) (5) 

The free field hr' satisfies a reduced Dirac equation 

-ivdh (v) (x) = (6) 

and gauge interactions are incorporated simply by replacing the partial 
with the covariant derivative 

d - D = d + [A,...] (7) 

For each value of the classical velocity v, there is an independent reduced 
field. The dynamics of are expressed by the effective Lagrangian 

S£ (v) = h (v \x) ivDh (v \x) (8) 

Which is known as the Heavy Quark Effective Theory (HQET). 

A lattice version of this theory requires an appropriate 
discretization of the reduced Dirac equation[4]. To implement on the 
lattice the property that, in the M -* °° limit, heavy quarks propagate only 
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forward in time drives one to an asymmetric discretization of the 
derivative operator. The most convenient choice is to use a symmetrical 
first difference for the spatial derivative but an asymmetric forward 
difference for the time derivative, a form often referred to as forward- 
time centered-space (FTCS). 

v [U(x,x+f) S(x+f,y) - S(x,y)] 

3 - iv (9) 
+ ^[U(x,x+p)S(x+p,y) - U(x,x-p)S(x-p,y)] = 5(x,y) 

This equation effectively defines the lattice heavy quark theory. 

The lattice reduced heavy quark propagator is the solution to 
this equation. On each lattice in a Monte Carlo ensemble it is obtained 
by direct forward recursion, which requires only a tiny amount of 
computation. 

If the lattice interaction strength is zero, [ U(x , x + p ) = 1 ] , the 
reduced propagator is proportional to 

e [-E(k)t + ik-x] (JO) 



where the energy E is given by 



E(k) = -In 



(11) 



This is the lattice approximation to the continuum behavior 



E(k) = lim [\/m 2 + {Mv+k) 2 -M\j\ +v 2 



v-k 



(12) 
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3. LATTICE HEAVY QUARK RENORMALIZATION 

This is surely not the place to review all the standard material 
on renormalization or even the renormalization of the heavy quark 
effective theory! The thorough review of Neubert covers the latter, 
including such subtleties as the matching of the effective heavy quark 
theory to QCD and the renormalization of composite operators, which 
enter into the calculation of weak decay rates[5]. 

On the lattice there is, however, a new renormalization which 
has no continuum analogue. It is the renormalization of the heavy quark 
classical velocity. Recall that in introducing the classical velocity, we 
described it in terms of the division of the heavy quark's momentum 
into a large fixed (classical) part and a dynamical part whose magnitude 
is set by the natural energy scale of QCD. We then used this division to 
arrive at a reduced Dirac equation for the heavy quark field. In principle, 
these two usages of the heavy quark classical velocity could logically be 
different. 

In the continuum, the values of the velocity that appear in the 
Dirac equation and as a parameter in the separation of the total four- 
momentum are always the same. To see this, it is sufficient to realize 
that v is the only four-vector parameter in the heavy quark theory. 
Therefore, the classical velocities appearing in those two contexts must 
be proportional to each other, and since they are by definition both 
normalized to 1, they are equal. 

However, on the lattice this is not the case[6]. A difference is 
possible because of the reduced symmetry of the lattice relative to 
continuum space-time. In the lattice heavy quark theory there are many 
linearly independent "four-vectors" that can be made from the 
components of a given four-vector v, where by a lattice four-vector is 
meant any quantity that has the same transformation properties as v 
under the lattice rotation-reflection group. The simplest examples are 

, 2n + 1 2n + 1 2n + 1 2n + 1 *. , 

the vectors (v ,v t ,v 2 ,v 3 ), whose components are the 
(2n+l) st power of the corresponding component of v. 
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Let us take as the definition of the bare, or input lattice classical 
velocity the quantity that appears in the reduced Dirac equation. In 
perturbation theory, this is the quantity that appears in the bare 
propagator. The physical classical velocity is the quantity that appears 
in the break-up of the four-momentum into a part proportional to the 
heavy quark mass and a part which is finite in the heavy quark limit. On 
the lattice, the components of this physical classical velocity could be 
proportional to any linear combination of odd powers of the components 
of the input classical velocity, subject only to an overall normalization 
condition. Since the possibility of a renormalization of the classical 
velocity arises purely because of the discretization of space-time, the 
actual shift is very sensitive to the details of the how the heavy quark 
theory is implemented on the lattice. 

Renormalization is most frequently discussed in perturbation 
theory, but in the case of the classical velocity, as we shall see, it is 
probably not at all reliable. To keep the context clear, when referring to 
the physical or renormalized classical velocity, I will use different 
symbols depending on whether it is determined by simulation or by 
lattice perturbation theory. Throughout this lecture I will use the symbol 
v (phys) f Qr ^ non perturbatively determined physical classical velocity, 
and the symbol v^ ren ' for the same quantity when it is calculated in 
perturbation theory. For the input or bare classical velocity I will use 
v (input) ^ Qr s j m piy Vj except when the perturbative context is to be 
emphasized, in which case the symbol v^ are ' will be used. 

In perturbation theory it is most natural to directly evaluate the 
residual heavy quark propagator, but in simulations is both more natural 
and more convenient to simulate the residual propagator of a composite 
particle containing a heavy quark. In order to simulate the heavy quark 
propagator, we would have to apply a global gauge-fixing procedure, 
and to exploit residual lattice symmetry we would have to choose a 
covariant gauge, such as the lattice Landau gauge. In perturbation 
theory this is all accomplished automatically by the choice of the gluon 
propagator, but non-perturbatively, gauge fixing would demand rather 
more computation than the calculation of the heavy quark propagator 
itself. 
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The physical classical velocity of a meson made of one heavy 
and one light quark is the same as that of its heavy quark component. 
The reason for this is because the mass difference between the heavy 
quark and the composite meson remains finite in the heavy quark limit, 
but corresponds to a finitely different breakup of the total 4-momentum 



Mv + p - (M + m)v' + p' 



(13) 



In the heavy quark limit, both v - v ' and p p' vanish like mIM, 
where m is the difference between the composite particle mass and the 
heavy quark mass M. 

Let us define the physical classical velocity, v^^, through the 
continuum physical residual energy of a particle containing a heavy 
quark, relative to the energy of its quark component. 



E (v \p) = lim y/(M + m) 2 + ((M + m)v <J,hys) + p) 2 - Mv C 



(phys) 



(phys) 



(14) 



Aphys)/ (phys\ 
i IV Q > 



From this relation we have 



,(phys) 



dE (v \p) 

d p 



(15) 



On an infinite lattice the momentum is continuous (though bounded), so 
this definition is quite sufficient. On a finite lattice we must specify a 
finite difference approximation to the momentum derivative. We will 
use the symmetrical first difference: 



(phys) 



A£ (v) (p) 



Yp~. 

* mm 



(16) 
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Here = 2n/Na is the smallest finite momentum representable on a 
lattice with spacing a and N sites to a side. 



4. RENORMALIZATION OF THE CLASSICAL VELOCITY 

The simulation of the classical velocity requires the computation 
of the function £^(p), which controls the asymptotic rate of falloff of 
a heavy particle's propagator M (v) (f,/f). For fixed v and p, the 
asymptotic behavior of the propagator is 

M (v) (t,p) ~ C (v \p)e~ EM( P )t (17) 



The energy E^(p) is the asymptotic coefficient of the Euclidean time t 
in the logarithmic derivative of M^' with respect to p, evaluated at zero 
momentum. 



dM (v >(t,p)/dp i |, 



p = 



M (v) (t,p =0) 

dC^\p)ld Pi \ pmQ 
C (v) (p=0) 

dC^(p)/d Pi \ ?m0 



C (v) (p=0) 



dE {v \p) 



d Pj 



p = 



v i 1 



(18) 



To adapt this discussion to the actual situation of finite-extent 
lattices, the only modification required is the replacement of the 
continuum momentum derivative by a lattice approximation, which we 
again take to be the symmetrical first difference on the finite Fourier 
transform lattice. 
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dM (v \t,p) 



dp, 

1 

2 prm 



\M {v \t,p) 



M 



p = 

M (v) (t,p = f Pima ) 



,7 = 

M (v) (t, p 



(19) 



fp ■ ) 

r min-' 



4.1 Expansion in the Bare Classical Velocity 

The structure of the calculation of the propagator M (v) (t,p) 
and the quantities derived from it is brought out by expanding in powers 
of the components of the bare classical velocity: 

M^(t,p) = £ ^^^M(t,p,m) (20) 

nij ,m 2 ,m 3 

The gauge invariant propagator satisfies lattice rotation- 
inversion symmetry conditions. For the case of a scalar meson, the 
propagator is invariant under simultaneous lattice transformations of v 
and p. To express this in terms of the coefficients M(t,p,m), recall that 
the 48-element three-dimensional lattice rotation-reflection group is 
equivalent to the group of permutations on the three axes times 
independent inversion groups for each axis[7]. The coefficients are then 
invariant under simultaneous permutations of the components of 
p = (pi, P 2, P j) and m = (rr\ , , ). In addition, M(t,p,m) is an 
even or an odd function of Pi depending on whether m ■ is an even or 
an odd integer, respectively. Note that because of the lattice symmetries, 
the derivative or symmetrical first difference of M (v \t,p) with respect 
to pi at p = must be odd in v { and even in the orthogonal 
components v. (j * i) . 

The basic elements of the simulation are the coefficients 
M(t,p,m) which appear in the expansion of in powers of the input 
classical velocity. Since the physical classical velocity is directly 
extracted from the logarithmic derivative of M^ v \t,p), we expand it in an 
analogous series in the input classical momentum v : 
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= £ '! V' 3 3 R (l \t,m) (21) 

M (v \t,p=0) 



Each of the expansion functions R (l \t,m) is asymptotically linear in t, 
and the negative of its slope is the corresponding term in the expansion 
for v- phys) in powers of ? : 

R (i) (t,m) ~ Z (i \m) - c (i \m)t 

^ = j2 *> 2 %r c «Km (22) 

m 1 ,m 2 ,m 3 



4.2 The Heavy Quark Propagator Coefficients 

A heavy meson propagator is computed in lattice simulations 
by calculating on each lattice of an ensemble a heavy quark propagator 
and a light quark propagator, and averaging the product over the 
ensemble. The light quark propagator is computed by solving the Dirac 
equation on the lattice. There are many methods in use, including the 
original Wilson formulation with a term proportional to the lattice 
spacing added to control the fermion doubling artifacts that, in one guise 
or another, must be overcome in all simulations. Since our focus here is 
those aspects of the simulation that are peculiar to heavy quarks, we will 
discuss the light quarks no further. 

The classical velocity enters through the equation for the 
reduced heavy quark propagator, Eq. (9). If we divide through the lattice 
Dirac equation by v , we see that each component of the three-velocity, 
v i = v/vq plays the role of a transverse hopping constant. That is, the 
reduced heavy quark propagator gets one factor of ?■ for each unit 
displacement in the + f direction. Since the propagator at any 
particular site on the n* time slice depends only on its value at the same 
spatial site or sites displaced by one transverse lattice spacing on the 
previous time slice, after n time steps the propagator is an (n - l) st 
order polynomial in the components f i . Since the classical velocity 
enters only through the heavy quark propagator, each coefficient in the 
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series for the meson propagator involves only one term in the heavy 
quark propagator series, that with just the same value of m. 

The calculation of the heavy quark propagator in simulations 
is greatly facilitated by exploiting this structure and expanding the 
propagator in a power series in ? ■ : 

S<%,*) = £ S(t AM) (23) 

m 1 ,m 2 ,m 3 

The site x is measured from the starting location. The computation of 
the coefficients in this polynomial is highly efficient. The value of the 
index m i is the maximum transverse lattice displacement of the heavy 
quark propagator in the i ih direction contributing to that coefficient. 
More precisely, the coefficient S(t,x,m) is non-zero only on the sites 

x i = -m i , -m ( . + 2 ,...,+ m i (i = 1,2,3) 

3 (24) 
h m t < t-\ 
i= i 

Since each transverse hop in the computation of S ( t , x , m ) can occur 
at any time between the initial and final times, the relative growth of 
these coefficients with t will have a factor 

S(t,x,m) - r. + ™ 2 + ™3 (25) 



4.3 Variational optimization of the meson creation operator 

It is well known that a gauge-invariant heavy-light meson field 
consisting of component fields at the same point in space-time is a poor 
choice in simulations, because many time steps are needed for the 
excited meson states to die off and for the ground state to dominate its 
propagator. While any smeared operator containing the reduced heavy 
quark field can be expected to improve the rate of convergence, in actual 
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simulations of quantities defined via composite fields, the limiting factor 
in the precision of the final result is usually how quickly the ground state 
in a given sector dominates the propagator. This is so because the 
statistical precision of propagators deteriorates rapidly with increasing 
lattice time, and so it is crucial that the propagator reach its asymptotic 
form in as few time steps as possible. 

Draper, McNeile, and Nenkov[8] have formulated this 
requirement variationally. Their idea is as follows. In any sector defined 
by a given set of quantum numbers and a residual momentum, the 
ground state contribution to a composite field propagator has the slowest 
rate of decay. Therefore, the requirement that the composite field be 
chosen so that its propagator reaches its asymptotic form in as few time 
steps as possible is equivalent to the requirement that the state created 
by the composite field have as large a component as possible along the 
ground state. This in turn is equivalent to the requirement that for any 
given time separation, the composite field is chosen so that the 
magnitude of its propagator is as large as possible. 

Let us consider composite operators that, in the Coulomb gauge, 
as functions of time and residual momentum have the form 



(26) 



where q(x) is a light quark field, i|j is a relative coordinate weighting 
function and the sum goes over sites y on the same time slice as x. For v 
= 0, t|j is effectively the wave function of a composite particle, but in 
general the weighting function does not have this interpretation. The 
propagator of this composite field is 

M^(t,p) -J2VS;(y)K^(t,y,y')^(y') (27) 



where K is the composite field propagator for fixed total three- 
momentum and arbitrary initial and final separation. 
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K$(t,y,y') = £ (s(x,x'=0)S (v \x + y,y')) (28) 

The average is over the ensemble of lattices. While it might seem 
counter-intuitive to take the location of the light quark as nominal 
coordinate of the meson field, it is by far handiest to do this. The 
requirement that the normalized propagator be maximal on a given time 
slice implies that t|/ v - ) is that eigenvector of K^ v ' with the largest 
eigenvalue: 

£ K$(t,y,y')^;](t,y') = k™p)y™(t,y) (29 ) 

i 

The largest eigenvalue is, in fact, the value of the optimal meson 
propagator on time slice t. 



The solution to this eigenvalue problem is facilitated by 
expanding the kernel, eigenvalue, and eigenvector in power series in the 
components of the classical velocity p. . 



K^\t,y,y') = £ 'i V»3 



i (v)/ n m, m 2 m 3 , ^ 

A ^(0 = L ff l V 2 *3 (31) 

( v ) m, m 2 m, _ 



m 1 ,m 2 ,m 3 



In terms of these expansion functions, the eigenvector/eigenvalue 
problem takes a familiar form. It has the structure of a perturbative 
treatment of an eigenvalue problem, with the components of the input 
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classical velocity playing the mathematical role of coupling strengths. 
In particular, it is only the static, zeroth order equation that requires 
finding the eigenvalues of a matrix, and it is only the zeroth order kernel 
matrix that must be inverted. Both the final time t and the residual 
momentum of the meson p are fixed parameters. 

The dependence of the zeroth order kernel on the residual 
momentum may be explicitly displayed, so that once the p = problem 
is solved, the results for non-zero p can be immediately read off. The 
zeroth order kernel 



K [p p,y,y,m = Q) = £ e ~ ,/7 ' A >(x,x' = 0) S(x +y,y,m = 0)) (32) 



is quite simple because the residual heavy quark propagator has no 
transverse hops. Thus the sum over x has only one contribution, at 
x = y' - y , which gives the momentum dependence as a diagonal 
similarity transformation of the zero-momentum kernel: 



K^(t,y,y') = e i ^ 7 ^(y-y,x' = 0)5 (v) (f,y';0,y')>^" , ' /7 ' f (33) 



The zeroth order eigenvalue is thus independent of the residual 
momentum, as it should be, and the residual momentum dependence of 
the zeroth order eigenvector is simply a phase. 



The non-trivial equation for the p = zeroth order eigenvalue and 
eigenvector is: 



-V 



ty m (t,y,m = 0) 
X m (t,m = 0) 



e i P^ m (t,y,m = 0) 
X m (t,m = 0) 



(34) 



Y / K [Q pJ,y',m = 0)\\i m (t,y',m = 0) 



f 



(35) 



= X m (t,m = 0) i|/ r01 (f,y,m = 0) 
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The higher terms in the expansion of the eigenvalue equation 
Eq. (29) are semi-negative-definite inhomogeneous equations. 
Suppressing the matrix labels (y , y ' ), the equation for the general term 
is 

[K m (t,m = 0) - k m (t,m = 0)]ty m (t,m.) 

(36) 

in ' 

The sum goes over all m' satisfying < m ' < m ., but excluding m' = 0. 
These are actually a hierarchy of inhomogeneous equations for the 
expansion coefficients of both the eigenvalue and the eigenfunction. The 
equation for the m th terms involves all the previous ones, that is those 
with m. < m p but subject to the strict inequality 2m' ( . < Sm r Note that 

even though i|j [0] (t,y,m = 0) is invariant under cubic transformations 

of y, all higher, m * terms are always asymmetric. 

All of the inhomogeneous equations for varying m ' have the 
same kernel, 

[K m (t,y,y',m = 0) - k m (t,m = 0) 8 ?>? ,] = 

(37) 

e^y[K m {t,y,y',m = 0) - X m (t,m = 0) b^]e 

The component of the equation for the general term, Eq. (36) along 
^[p\ ( f ' = ) gives the eigenvalue coefficient: 

in 

^ m (t,m) = Y, ^ m {t,m = Q)K m {t,m')^ m {t,m-m') (38) 

m ' f 

The kernel restricted to the orthogonal subspace is negative definite and 
so can be inverted nonsingularly on that space: 
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y m (t,m) =-{Q m [K m (t,m = 0) - k m (t,m = 0)]Q m y 1 

m (39) 

x £ l K m (t> m i') - k m (t>m.')]y m <.t,m.-m.') 

m' f 6 



where 

Q m ( y,f ) = « v- - ^roi ( * & * = ) i|4 ( f ,f, m = 0) 



(40) 



is the projector on the space orthogonal to the ground state. The matrix 
inverse is nonsingular on the image space of j. Note that only a 
single matrix needs to be inverted throughout this process. 

The successive terms in the expansion of the weighting function 
each have different large Euclidean time behaviors. Specifically, the 
relative growth with t of the terms in the expansion of the heavy quark 
propagator with different values of m (Eq. (25)) results in the same 
relative t dependence of the terms in the expansion of the kernel (Eq. 
(28)). Since tis a fixed parameter in the whole hierarchy of expressions 
for the wave function components (Eq. (39)), the coefficient 
[p] ( t , y , in) grows with t like 

y m (t,y,m) ~ V 1 2 3> (41) 



5. SOME SIMULATION RESULTS 

This is as far as we can go analytically. To actually solve the 
variational problem, we must specify over what class of functions we 
will maximize the normalized propagator Eq. (27), and then numerically 
solve the resulting matrix eigenvalue problem. 

A simulation was carries out on an ensemble of lattices and 
Wilson light quark propagators created by the Fermilab ACP-MAPS 
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Collaboration^]. The ensemble consisted of 42 lattices of size 24 3 x 48 
with lattice coupling P = 6.1 along with Wilson quark propagators with 
hopping constant k = .154. Heavy quarks propagate in only one 
Euclidean time direction, but which direction is conventional. Each 
lattice therefore provides effectively two independent configurations, one 
each for association of the direction of time propagation of the heavy 
quark with the nominal forward or backward direction of lattice time. 

5.1 The Weighting Function 

As the space of functions over which we maximize the 
composite propagator, we choose those which are non-zero only within 
a box with (2N W + 1) sites on a side, and invariant under lattice cubic 
rotations and reflections within it. We vary N w to find how sensitive 
results are to the size of the function space. As a consequence of the 
symmetry condition, the rank of the square matrix to be inverted is 
(N w + l)(N w + 2)(N W + 3)/6, which is fortunately much less than (2N W 
+ l) 3 , the value it would take without imposing the requirement that the 
ground state of the composite particle be a lattice scalar. 

Figure 1 shows how the quality of the solution improves with 
increasing N w . What is plotted is the effective mass of the ground state 
as a function of time slice, for box sizes ranging from N w = to 
N w = 6. The effective mass means, in this context, the log of the ratio 
of the composite propagator on two successive time slices. The 
statistical errors on each point are determined independently for each 
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Figure 1 — Convergence to an effective mass plateau of a static 
heavy-light meson propagator optimized within domains containing 

(2N W + l) 3 sites. 



value of N w and each pair of adjacent time slices using a single lattice 
ensemble of lattices. The curves are from two-exponential fits to the 
underlying propagators. 

There are many lessons to be learned from these plots, which 
nicely show the trade-offs between finite statistics, which favors short 
time differences, and the exponential damping of excited state 
contributions to propagators, which favors long separations, as well as 
between computing resources, which favors small boxes, and rapid 
convergence, which favors large ones. For example, it is well known 
that asymptotically in lattice time only the ground state with given 
quantum numbers contributes to a propagator with those quantum 
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numbers. This is the source of the statement that any operator with the 
correct quantum numbers is an acceptable interpolating field. However, 
it is clear from the two exponential fits that the asymptotic mass seems 
to depend quite strongly on the size of the box. The reason for this 
apparent contradiction is that, because of the deterioration of statistics 
with increasing lattice time, the weighting of successive points in a fit 
falls off with lattice time. Effectively, only a few, small-separation 
points actually determine the fit and its attendant asymptotic mass. 
Since the contribution of the high mass states to the composite 
propagator at any given lattice time is less for larger boxes, the mass 
derived from each fit is different. The values derived from each fit 
converge to the true value as the box size increases. Using larger 
ensembles of lattices improves statistics, which ameliorates this effect, 
but can never eliminate it entirely. 

With the above remarks as background, let us note several 
characteristics of the results. The first thing to note is that the local 
(N w = 0) operator never reaches a plateau. That is, without some sort 
of smearing procedure for enhancing the convergence of the composite 
propagator, this program for non-perturbatively computing the 
renormalization of the classical velocity would fail. Another observation 
is that the statistical errors, shown by the error bars, are not a true 
measure of the uncertainty in the values of the effective mass. To see 
this requires looking closely at the graph. Notice that there are parallel 
undulations of the data points for the different sized domains relative to 
the fits. It is the magnitude of these undulations, not the statistics of each 
data point, that is a meaningful measure of the precision with which the 
effective mass is determined. The reason that the undulations move in 
parallel is that they are reflecting correlations from one time slice to the 
next within the ensemble of lattices. To within that precision, no evident 
plateau is found for N w = or for N w = 1 (27 sites). For N w = 2 
(125 sites), a plateau is reached by about t = 4. For the larger domains, 
N w > 2 (7 3 to 1 3 3 sites), a plateau is evident by time slice t = 3. In 
fact, given the limited precision indicated by the parallel undulations, 
one cannot be sure that a plateau is not reached even more quickly. 
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5.2 The Classical Velocity 

Let us now see how these ideas work in calculating the shift in 
the classical velocity. Because this shift is a finite renormalization 
which vanishes in the continuum, as opposed to the familiar situation in 
which renormalizations result from continuum divergences, it is very 
interesting to see whether the shift is a large effect or simply a matter of 
principle without much practical import. 

The starting point for the optimization of the composite meson 
operator described in the previous section is the propagation kernel 
defined by Eq. (28), or more precisely the coefficients K^(t,y,y',m) of 
its expansion in powers of v. These are the quantities that are the 
immediate output of a Monte Carlo simulation. We evaluate these 
coefficients for relative coordinates in a box of with 5 sites on a side 
(N w = 2). This is a compromise value. The choice fy = 3 would 
yield a more rapid approach to a plateau, but would require a great deal 
more computer memory. Higher values of N w would be prohibitive in 
their use of computer resources. 

The order of battle is as follows. We solve the zeroth order 
eigenvalue equation, Eq. (35) for its highest eigenvalue and its 
associated eigenvector with p = 0. We then use Eq. (34) to obtain the 
zeroth order eigenvector for p = (1,0,0) and its permutations, and then 
iterate Eqs. (38) and (39). From the classical velocity expansion of 
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Figure 2 — The zeroth and first order components of the optimized 
heavy-light meson propagator, showing degree of dominance by the 
ground state (The right axis is the fractional saturation.) 



Eq. (30), or equivalently of Eq. (27), we find the leading term in the 
classical velocity expansion of the variationally optimized heavy-light 
meson propagator. 

The zeroth and first order terms in the expansion of the heavy- 
light meson propagator A/^ that enter into the determination of the 
physical classical velocity are shown in Figure 2. Since terms odd in m i 
are odd in p { , the values given for M(t,p = (\..),m = (odd..)) are 
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exactly the same as = 2n/Na times the lattice approximation to the 
momentum derivative AM/Ap ; .|, = of Eq. (19). 

There is a lot of information in this graph, and it is worthwhile 
examining it in a bit of detail. Note the high precision of the leading 
terms and their convergence to the expected asymptotic forms. The error 
bars on the computed propagator values are from the single-elimination 
jackknife analysis, and the error bars on the percentage saturation are the 
fractional errors on the static propagator. The curve through the static 
propagator points is a double exponential fit: 

M(t,p =0,m = 0) ~ Cl e~" ht + c 2 e~" ht (42) 

The masses in the fit are m l = .580 + .012 and = .961 + .031, and 
the ratio of the coefficients is c 2 /c j = 2.20 ± . 1 3 . The higher mass term 
effectively models all of the higher mass effects in the static propagator. 
The straight line plotted under the static propagator is the first term from 
the above fit. The asymptotic mass is, of course, very slightly lower 
than that obtained from the time-slice to time-slice fall-off the 
propagator below t = 10. The computed saturation values are the single 
exponential component of the two-exponential fit divided by the 
propagator values, and the saturation curve is the ratio of the single 
exponential to that fit. 

The curve through the computed values of the first order 
component of the meson propagator is a fit to the expected asymptotic 
form of this term, a linear function of lattice time times the two- 
exponential fit to the static propagator, 

M(t, p = (1,0,0), m = (1,0,0)) ~ (a+bt)( Cl e~ mit + c 2 e~" hJ ) (43) 

The fits uses the values of the coefficients c ■ and m ■ determined from the 
static propagator. As we noted before, the asymptotic slopes of 
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Figure 3 — Leading tj Expansion Function in the Logarithmic 
Derivative of the Meson Propagator 



the functions R (l) {t,m) of Eq. (21) are the negatives of the coefficients 
of the expansion of the physical classical velocity in powers of the input 
classical velocity. 

Figure 3 is the graph of the negative of the leading function, 
R (l \t, m . = 1 , m. , .= 0) . The indicated statistical errors are obtained from 
a single-elimination jackknife procedure. This leading term in the 
logarithmic derivative of the meson propagator apparently achieves its 
"asymptotic" slope on the very first time slice on which it is non-zero, 
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t = 2. The statistical error in the fitted slope, .885, is + .017, and so to 
leading order the shift is v ( phys) / v ( mput) = .885 + .017. Corrections 
come from the cubic and higher powers of v- mput) m the expansion of 

(phys) 

It is both interesting and disturbing to compare this result to a 
perturbative, one-loop calculation of the renormalization of the classical 
velocity [6]. For the lattices used in the simulation, the coupling constant 
is g 2 = 6/p = 6/6.1 , and the linear term in the shift of v i is -.23345569. 
The comparison of the leading terms is 



The disagreement between the two methods is quite substantial, because 
the perturbative calculation evaluates the difference of the above ratio 
froml. There are, of course, higher loop terms in perturbation theory, 
but to bring these two results into concord, higher loops would need to 
cancel just about half of the leading term. 

This is by no means the first example of a serious discrepancy 
between lattice simulation and lattice perturbation theory. However, the 
rescaling of the classical velocity is different from other 
lattice/continuum matching processes is that it seems to be exacerbated 
rather than improved by the "Tadpole Improvement" procedure of 
Lepage and Mackenzie[10]. 

A principal quantitative effect of Lepage and Mackenzie's 
analysis to introduce a scaling between the couplings used in simulations 
with those used in the perturbative analysis. It arises because on the 
lattice there are tadpole corrections to the gauge field, and to cancel these 
lattice artifacts, one replaces each link variable U^n) by that link 
divided by its average value, U (n)/u , where the average value of the 
link is taken to be the gauge-invariant expression 



.885 + .017 simulation 



(phys) , (input) 
i < V i 



.76654431 one loop 



(44) 



-{TrU) 
3 



1/4 



U. 







(45) 
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and □ denotes the product of the link variables around an elementary 
lxl plaquette. For the Wilson action used to generate the Fermilab 
lattices, this effects the change 




and so P is identified with 6/g 2 u Q 4 rather than with 6/g 2 . Since u Q is 
necessarily less than one, this always has the effect of increasing the 
value of the effective perturbative coupling corresponding to a lattice 
simulation at a given value of p, thereby increasing the discrepancy of 
the perturbative estimate of the classical velocity shift from the result of 
the simulation. Perhaps the fact that tadpole "improvement" worsens 
the disagreement between the two calculations is related to the fact that 
the classical velocity renormalization has no divergent part and vanishes 
in the continuum, but this is just conjecture. 

6. THE FUTURE - THE HEAVY QUARK FORM FACTOR 

The paradigmatic heavy quark process is the weak decay of one 
meson containing a heavy quark into another, such as the B meson 
decays 

B - Dlv , D*lv 

(47) 

(Dv'| cy m (1 ~iy 5 )b\Bv) ~ l(vv') 

All the dynamical information is contained in the function E, , known as 
the Isgur-Wise function, or the heavy quark universal form factor. It is 
a function only of the dot product of the classical velocities, and is 
normalized to unity for forward decay (£(1) =1) because the same form 
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Figure 4 — A heavy meson decay in which one heavy quark emits a 
weak current and transforms into a different heavy quark 

factor appears in the B meson matrix elements of the conserved vector 
current. The amplitude is contained in the reduced three-point function 

G (v *' VD) (x ,z ,y ) = (rr S (y,x)r D 5 (VD) (x,z)r y 5 (V (z,y)r B ) (48) 

which is schematically illustrated in Figure 4. 

The Isgur-Wise function may be extracted from the three-point 
correlation by dividing out the normalized single meson propagators. A 
better way is to take advantage of the normalization in the forward 
direction and use the three-point function at v D = v B as normalization, 
giving 
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£(vv') | 2 = lim 



•o 



) G«- v \x Q ,y Q ,z t 



•o 



) 



(49) 



G (v ' v) (x ,y ,z 



) G( v '- V '>(x ,y ,z, 



■o 



) 



This form has the advantage that the overall normalization of the three- 
point function drops out completely. Only that part of the 
renormalization that is classical-velocity-dependent need be included. 
Unfortunately, this renormalization is not at all simple to compute by 
simulation, but it is crucial to deriving physical results from this 
analysis. The reason is that the slope of the Isgur-Wise function at the 
kinematic origin is needed in phenomenological applications. This is a 
key next step in the study of the lattice heavy quark theory. 

The renormalization of the classical velocity affects the 
calculation of the slope of the Isgur-Wise function using the lattice 
HQET simply because it changes the heavy quark momentum scale. Of 
course, only the linear, multiplicative shift enters into the slope of the 
Isgur-Wise form factor. The cubic and higher order terms that we have 
not calculated only affect its higher derivatives. The precise effect is 
seen by expressing the standard kinematic variable w = vv' - 1 in terms 
of the bounded classical velocity v 

w = vv' - 1 = i(» - v') 2 + 0((v,v') 4 ) (50) 

We see that it is second order in v. Thus the rescaling v - has 
the effect of increasing the computed slope of the form factor with 
respect to vv' at the origin by the square of the first order rescaling 
coefficient, 



1 



1.28 ± .05 



C (,) (fflfl) 2 



(51) 
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7. CONCLUDING REMARKS 



Let us sum up what we haave learned about the heavy quark 
effective theory on the lattice. We showed how to formulate the Isgur- 
Wise limit in Euclidean space, and examined some of the properties of 
the a lattice discretization of the reduced Dirac equation for the residual 
heavy quark propagator. We examined in detail the renormalization of 
the classical velocity in this treatment. The origin of this 
renormalization, which is purely a lattice effect and is not present in the 
continuum, is the reduction of Lorentz invariance to hypercubic 
symmetry on the lattice. We reviewed the result of a lattice simulation 
of the first-order, multiplicative shift, and obtained 

v (phys) = (_ 885± _0 17 ) f . + 0(vf,V.vh 

I I I I J 



We discussed how the lattice HQET provides a method for calculating 
the nonperturbative functions that arise in the heavy quark limit of QCD, 
especially the Isgur-Wise form factor, but noted that one must first carry 
out a simulation of the renormalization of the composite weak current 
operator. Finally, we noted that the renormalization of the classical 
velocity will have a noticeable effect in such calculations by effectively 
shifting the scale of the Isgur-Wise kinematic variable. 
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